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Here, we develop a theory of radiative heat transfer based on an equivalent electrical network 
(N- representation for the hot material slabs in an arbitrary multilayered environment with arbitrary 

distribution of temperatures and electromagnetic properties among the layers. Our approach is fully 
equivalent to the known theories operating with the fluctuating current density, while being signifi- 
cantly simpler in analysis and applications. A practical example of the near-infrared heat transfer 
^ ' through the micron gap filled with an indefinite metamaterial is considered using the suggested 

O ■ method. The giant enhancement of the transferred heat compared to the case of the empty gap is 

shown. 
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I. INTRODUCTION 



There are two most important heat exchange processes known: thermal conductivity, associated with collective 
oscillations of atoms (i.e., phonons) or electrons in solids including metals, or with convection in fluids or gases, and 
radiative heat transfer which is associated to the electromagnetic radiation produced by thermally agitated atoms, 
e.g., the black body radiation. In this work we concentrate on the latter process which can be dominant when the 
bodies that exchange heat are separated by gaps (empty or filled with a heterogeneous material that weakly conducts 
the heat). 

^ '■ As is known, thermal radiation results from fluctuations of charge and current density in matter. This phenomenon 
O . is governed by the fluctuation-dissipation theoremi. (FDT) that relates the mean-square fluctuations of a physical 
'■ quantity to the dissipation associated with the dynamics of the same quantity. Because in dielectrics the loss is 
I represented by the imaginary part of the permittivity, which may be in turn related to the effective conductivity of 



the material, the FDT requires the volumetric current density within a dielectric or conducting body to fluctuate. 

■ The FDT constitutes the basis of the present day radiative heat transfer theories which deal with the fluctuating 
\^ , currents and treat them as the principal source of thermal radiation. Such a picture places the radiative heat transfer 
lO ■ calculations in the framework of classical electromagnetic theory based on the macroscopic Maxwell equations. 

^^O , Indeed, the well-known theory of radiative heat transfer through narrow vacuum gaps by Polder and van Hove^ 
belongs to this class. Similar techniques have been recently developed for the cases when radiative heat transfer 
through micron and even submicron gaps is assisted with nanostructured metamaterials (see, e.g., in Refs. In 
CN ■ these approaches the fluctuating current density in the two neighboring media is first obtained from the FDT, following 
the methodology introduced by Rytovi^ Next, the electromagnetic field produced by the fluctuating currents is 
calculated, after which the mean value of the power flow (Poynting vector) across the gap (fiUed with the metamaterial) 
is found either with the multiple reflection method,-'-'^ or with a more general transfer matrix approach?^ 

Historically, however, thermal agitation of fluctuating currents was first discovered by Johnson^ in electric circuits 

■ and the theory of such thermal fluctuations was developed by Nyquist^ using ideas which were natural for an engineer 
dealing with networks of lumped elements and transmission lines. Although this theory was a precursor to the FDT, 
it is still in wide use in the theory of thermal noise at radio and microwave frequencies. The famous Nyquist's result 
states that, in any linear passive two-pole (i.e., a single port device with an input represented by two electric contacts) 
operating at a temperature T the electric thermal fluctuations (in other words, the thermal noise) concentrated within 
a narrow frequency interval Av can be equivalently represented by the fluctuating electromotive force (EMF) e{t), 
with the mean-square of fluctuations 

= 4:e{iy,T)R{iy)Av, (1) 

where Q{v,T) — hv / [exp{hv / k^T) — 1] is Planck's mean energy of a harmonic oscillator with h and fce being the 
Planck and Boltzmann constants, respectively, and R{v) is the input resistance (real part of the input impedance) of 
the two-pole. The equivalent EMF is then understood as connected in series with the two-pole, which can be now 
considered noiseless. 
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The beauty of this result is in that no knowledge of the internal structure of the electric network is required, and 
that all the information is contained within just a single parameter: the real part of the frequency-dependent input 
impedance. In the terminology of FDT, the equivalent EMF in Nyquist's formula has the role of a generalized force, 
and the two-pole input impedance is related to the generalized susceptibility of the system. In other words, the 
Nyquist formula can be obtained by a direct application of FDT to the electric circuit.™ Note that the fluctuating 
current appears in this model only as a reaction to a finite number of lumped voltage sources. 

In contrast, in the thermal transfer theory based on the full- wave electromagnetic formulation through the fluctu- 
ating current density understood as a volume-distributed source, the internal structure of the interacting bodies has 
to be considered during most of the calculations. Thus, a thermal transfer problem appears in this formulation as a 
problem with infinite number of degrees of freedom. The final result, however, happens to be expressed in quantities 
that abstract away the internal structure, such as the reflection coefficients of material half-spaces in the theory of 
Polder and van HovCfS This observation suggests that introducing the distributed fluctuating current can be avoided 
in many cases. For example, in this work we prove that in stratified media the input impedance concept and the 
original Nyquist theory can be generalized and used not only in problems related to electric networks, but also in 
full-wave radiative heat transfer problems. This allows for a significant reduction in complexity of the analysis and 
makes the theory readily available for practical calculations. 

It has to be mentioned that an expression for the equivalent lumped EMF of the thermal noise for the case of a 
lossy material half-space was first derived by Rytov^ using an approach based on the fluctuating current density. This 
result constitutes the basis of the thermal noise theory of aperture antennas. ^'^"^ Rytov also worked on equivalent 
four-pole network representation of a hot material slabi^ Why then these results have not been widely used in the 
heat transfer problems? Perhaps, it is because the concepts of input impedance and the equivalent circuit description 
for full-wave electromagnetic problems are largely unknown among theorists working in the field of radiative heat 
transfer. In applied electromagnetics, however, it is well-known that stratified media can be very efficiently treated 
within the so-called vector transmission line theorj*ii (VTLT) which, in essence, assigns an equivalent transmission 
line network to every electromagnetic mode (propagating or evanescent) in the system. We would like to stress here 
that the VTLT is not an approximation: it is a direct consequence of the Maxwell equations when modal expansion 
is applied to the electromagnetic field in layered structures. The VTLT allows also for a systematic treatment of 
uniaxial and bi-anisotropic media. 

In this work we extend the VTLT in order to include the effect of the fluctuating current density within the layers. 
In contrast to a few numerical discrete-element models currently available from the literature (e.g., Ref. 11^ . the 
theory that we develop here is fully analytical. The generalized theory allows us to prove a complete equivalence 
between a volumetric multilayered structure and its circuit theory counterpart, which may be visualized as a chain of 
transmission line segments with equivalent fluctuating voltage sources connected at the ports. When concerned with 
the radiative heat transfer between the layers, we show that this equivalent network may be reduced to just a series 
connection of a number of voltage sources representing the fluctuating EMFs and equivalent impedances (each can be 
under different temperature), thus, recovering in this way the famous Nyquist result, generalized here to the full- wave 
electromagnetic processes in stratified media. Therefore, the calculation of the radiative heat transfer between the 
layers reduces in our theory to a number of equivalent circuit theory calculations, which are relatively simple and very 
similar to what is typically done when considering thermal noise in practical electric networksJ^^— 



II. IMPEDANCE REPRESENTATION OF THE CLASSICAL HEAT TRANSFER FORMULA 

In this section we establish a connection between the classic Polder-van Hove theory of radiative thermal transfer 
and its representation in terms of the input impedances of the material half-spaces. Formulas by Polder and van Hove 
for the density of the radiative heat flux (i.e., power flux of the thermal radiation) across the gap between two thick 
dielectric slabs [the geometry is deflned in Fig. [TJa); the slabs are approximated by half-spaces] read as (our notations 
correspond to Refs. [lO - llSh : 

oo 

/p- 
dw[e(c^,Ti)-e(c^,T3)]M, e(c^,T,) = ^^ — . (2) 

Here h = /i/(27r), and Ti is the absolute temperature of the i-th medium (i = 1, 3; i = 2 represents the gap). If Ti > T3 
the total radiative flux is directed from medium 1 to medium 3 and can we written as S't = S'l-j-s — 5*3^1 where 5*1^3 
is the heat flux produced by medium 1 and absorbed in medium 3 and S'3_i.i is the flux produced by medium 3 
and absorbed in medium 1. M from ([2]) is called the radiative heat transfer function^iSi It depends on the optical 
properties of all three media and may be written as the sum M = Mp + M^, where Mp and Me are contributions of 
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FIG. 1: (Color online) (a) - Illustration to the general problem formulation, (b) - A possible implementation of medium 
2 suggested in Ref. transforms the gap into a layered structure filled with so-called hyperbolic metamaterial (see, e.g., in 
Refs. i19h22) formed by carbon nanotubes or metal nanowires. Interdigital arrangement of nanotubes (nanowires) helps keep 
radiative heat transfer dominating over the thermal conductance through the gap. 



propagating and evanescent waves, respectively [note that in this paper we assume the time dependence exp(+jujt) 
with j = \/— 1] : 



Mp^^ I Np{cj,q)qdq, Np{u;,q) 



(i-|ri2(g,c.)p)(i-|r32(g,c.)p) 

4|I-e-2.7-/3<^ri2(g,w)r32(g,a;)P ' 



(3) 



M ^ [ AT( \ rl AJ ( \ Im[ri2(g,c.)]Im[r32(g,c.)]e-^l^l'^ 
Me^— Ne{u,q)qdq, Ne{uj,q) = — _2\Bidr ( vF~? ^I2~' 

ko 

where and r32 are the reflection coefficients of a plane wave harmonic having the spatial frequency (transverse wave 
number) q = |kt| and being incident from a lossless medium 2 (here, medium 2 is free space) to the surfaces of media 
I and 3, respectively, /3 — ^k^ — is the normal component of the wave vector in medium 2, and fco = '^^SoMo is 
the wave number in medium 2. Function ^'^(a;, q) is called the spatial spectrum of the radiative heat transfer function. 
This function was studied in Ref. [l^ for the case when photon tunnelling through the vacuum gap was enhanced 
by surface plasmon-polaritons excited at the gap boundaries. It was shown that the absolute maximum of N{uj,q) 
achievable at certain values of oj and q equals 1/4 (we discuss this limit with more detail in Section |VI|. 

Let us express reflection coefficients ri2 and r32 in ([3]) and (jl]) through wave impedances ^1.2,3 of the media in 
regions 1, 2 and 3. The wave impedance defines the ratio between the transverse electric and magnetic fields in a 
plane wave of a given polarization, i.e., in a given harmonic of the spatial spectrum (see, e.g., in Ref. ^ll'). In terms 
of the wave impedances, 

— ^2 

ri2 = ry , ry ^ * = 1> 3. (5) 

For isotropic dielectrics the wave impedances are given by the following expressions (see, e.g., in Ref. TT) for TM- waves 
and TE-waves, respectively: 

/7TM „ A ^TE ^0 



where rjo is the characteristic impedance of free space r/o — \/ l^o/£o, and Pi denotes the normal component of the 
wave vector in the i-th medium: /3i = \Jk'^ei — (p . If medium 2 is free space (as it is assumed in the classical theory 
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of Ref. '2') Ei = 1. For propagating waves q < and Z2 is real, for evanescent waves q > and Z2 is imaginary. It is 
important that the impedance representation ([5]) of the reflection coefficients is general for all spatial frequencies q. 
Substituting Eq. ([5]) into Eq. ([3]) we obtain after a rather simple algebra: 

N = IQRlRiR^ 

^ 4|(Zi + Z2){Z2 + ^3) + (^1 - Z2){Z2 - Z3)e-2,'3.d|2- ^ ) 

Here and below we denote Ri = Re(Zi) and Xi = Im(Zi). Substituting Eq. (O into Eq. ^ we easily deduce: 

^ ^ ^XlR,R,e-^\P^\^ 

" |(Zi + Z2){Z2 + ^3) + {Z, - Z2){Z2 - Z3)e-2|&l'i|2 ' ^ ^ 

In this case /32 has imaginary value and it is taken into account that Z2 — 3X2- Since X2 = for propagating waves 
and i?2 = for evanescent waves, Eqs. ([7]) and ([8|) can be unified into an expression suitable for both regions q < ko 
and q > ko: 

^ |(Zi+Z2)(Z2+Z3) + (Zi-Z2)(Z2-Z3)e-2j/5^'^|2- 

The heat flux density originated from the medium 1 and absorbed in the medium 3 can now be represented as 
integral over the complete spatial spectrum (both propagating and evanescent): 

00 00 

^1^3 = ^ Jdiojqdq e(w, Ti)N{lu, q). (10) 


As we show in the next section, Eq. © for the spatial spectrum of the heat transfer function iV(a;, q) can be derived 
from the VTLT in a way that allows for a straightforward generalization of the results of Polder and van Hove - to the 
case of stratified uniaxial magneto-dielectric media without a need to introduce the distributed fluctuating currents. 



III. RADIATIVE HEAT TRANSFER RESULTING FROM AN EQUIVALENT CIRCUIT APPROACH 

The possibility to express the radiative heat transfer through the wave impedances of the material layers provides 
us with an evidence that an equivalent circuit model can be formulated for this problem. Such a model is derived 
rigorously in Appendix directly from the Maxwell equations, which results in the VTLT generalized to include the 
effect of thermal fluctuations. Throughout this section, however, we will use simple physical reasoning when possible, 
in order to keep the derivations easy to grasp. 

As in the previous section, we decompose the fluctuating electromagnetic field into plane waves characterized with 
a certain polarization state, angular frequency uj, and the transverse wave vector kj. For generality, let us assume 
that all materials taking part in the heat transfer are optically uniaxial magneto-dielectric media. It is known that 
in uniaxial magneto-dielectrics the independent polarization states correspond to the TE and TM plane waves. The 
same holds for a multilayered structure composed of uniaxial magneto-dielectric layers (some of them can be isotropic 
or even vacuum gaps) under the condition that the anisotropy axes of all layers are aligned. In such a structure the 
two polarizations are completely independent and can be considered separately. In the following we assume that the 
anisotropy axis of the layers coincides with the z-axis which is perpendicular to the layers. 

In Appendix we prove that, in a given layer of the considered multilayered structure {i-th layer) being under the 
temperature Ti the transverse components of the time-harmonic fiuctuating electric and magnetic fields at the layer 
interfaces (labeled here with subscripts 1 and 2) are related as follows 




nixHl,\ /El,\ 1 



(11) 



where Z„„ — Z^^{ui, kt) are the dyadic Z-parameters of the chosen layer, and 2 are the vectorial fluctuating EMFs 
equivalently representing the thermal-electromagnetic fluctuations within the same layer, and ni_2 are the external 
unit normals at the interfaces of the layer. The meaninig of the factor Xj^fK^ is explained further in the text. Eq. pT|) 
generalizes the known result of the VTLT to the case of non-vanishing thermal fluctuations. When the right-hand 
side of (jlip vanishes this equation represents the definition of the impedance matrix of a passive material layer. 
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The result (|TT|) is obtained in a dyadic form and, thus, is apphcable to any polarization state of the electromagnetic 
field in an anisotropic (not only uniaxial) layer. However, when the states split into the TE and TM waves it is 
more convenient to work with scalar Z-parameters — elements of the 2x2 impedance matrix — which are defined 
separately for each polarization. For a slab of a uniaxial magneto-dielectric characterized by the permittivity dyadic 

e = e^It + ej'zoZo and the permeability dyadic JI — fij~It + mI'^oZq (here we understand these parameters as relative 

to the vacuum permittivity Eq and the permeability /j,Oj respectively, with /( being the unity dyadic in the transverse 
plane), the Z-parameters are (see, e.g., Ref . ITll) : 



^TE.TM 



^12 — ^21 — j- 



oTE.TM 



z. 



TE.TM 



^22 



tan(/3/ 



,TE,TM 



(12) 



where fSj^'^^^di is the electric thickness of the i-th layer, and the wave impedances of a spatial harmonic with wave 
vector k — {k^, ky, /3™''^^) are 



7TE _ ^ fcoMi'" 
— Vo- 



0, 



TE 



7™ - -n 

— Vo 



TM 



(13) 



where the propagation constants for the two polarizations are expressed through q 



k"^ as (see, e.g.. Illf): 



3TE 



\ 



PI 



TM 



\ 



g2) 



(14) 



Thus, when considering plane waves of fixed polarization and fixed transverse wave number q= |kt| such a slab is 
described by a 2 x 2 matrix of scalar Z-parameters, much like a four-pole network in the circuit theoryi^ To make the 
analogy complete, we may introduce the effective "currents" flowing into this four-pole network and relate them to 
the magnetic fields at the two interfaces of the slab as 



,TE rj- zq X kt , Y -tri ^ 



Iktl ^ ^^'^ lktl 

and the effective "voltages" at the input and the output interfaces of the slab 

.j,TM rr- kt j -,rt,TE _ rr- zq x kt . 



(15) 



(16) 



The factor \/Aq where Aq is the unit area in the transverse plane ensures that the complex power 
-Aq ni,2 ■ {^\a ^ -1-11,2 

polarization (TE or TM) relation ^ assumes the form 



TE.TM ^i.TE.TM" 
1,2 



TE,TM ^ jji,TE,TM \ trivially related to the complex Poynting vector of a mode. Then, for a given 



^11 ^12 

^21 ^22 



(17) 



The equivalent circuit that corresponds to this relation is shown in Fig. [51 The two EMFs at the input and the 
output of this circuit represent the effect of thermal fluctuations inside the chosen layer. The mean square amplitude 
of these equivalent sources is derived in Appendix using the approach of the distributed fluctuating current. However, 
one can apply the Nyquist theory directly to the electric circuit in Fig. [5] and obtain the same result. Namely, by 
disconnecting the load from the output of the four-pole network, i.e., setting = 0, we eliminate the contribution of 
e\ and obtain a simple two-pole network with the input impedance Z-^ — VI /I\ — Z\i. Therefore, the mean square 
amplitude's of the fluctuating EMF e\ within an angular frequency interval Aw is 

(4)^=2e(w,rORe(Zji)— (18) 
Repeating the same procedure while interchanging the roles of the input and the output one obtains that 



(ey2 = 2e(w,r,)Re(Z^2 



^ Aw 

TT 



(19) 
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FIG. 2: Equivalent four-pole network of a material layer under temperature T = Ti. 
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FIG. 3: (Golor online) An illustration to the calculation of the radiative thermal flux through a selected boundary z = Zi in a. 
multilayered structure. 



It is evident that, in general, the fluctuating sources e\ and 63 must be partially correlated, because they both 
represent the fluctuations within the same layer. Therefore, in calculations involving expressions which are quadratic 
in voltage and (or) current (e.g., power) one may also need the correlation function of these EMFs: {e\e\). As is 
shown in Appendix, this correlation can be presented as: 

Auj 



ie{el) = 2e(w,rORe(Zj2) . (20) 



TT 



The set of relations ([T5)) -(PU1) can be also obtained directly from the FDT. Indeed, if one identifies the charges 
9i,2 = ^1,2/ {3^) the state variables of the circuit depicted in Fig. [2] and the EMFs ei,2 as the random forces 
associated with the fluctuations, then the FDT demands that for the fluctuations concentrated within a narrow 
frequency interval Auj 



Z ZKB-L TT 

where a„m are the generalized susceptibilities such that qm = J2n Q:mn(w)e„. It is readily seen that {a~^)mn = j^Zmn- 
Substituting this into (|2T|) while taking into account the symmetry properties of Z„in we obtain (|18|) -([20 |) after 
dropping the irrelevant contribution resulting from the quantum zero-point fluctuations. 

The relations P^ - ipH)) together with (|17p written for both polarizations fully describe the fluctuations within a 
material layer. A structure formed by many layers can now be equivalently represented by a chain connection of many 
four-pole networks each representing a layer. Let us now select an arbitrary boundary between a pair of layers in a 
multilayered structure and find the radiative power flux per unit area of this boundary (Fig. [3]) . 

Because the fluctuating EMFs belonging to separate layers are uncorrelated, we may flrst consider only the sources 
which are located at z < Zi, where Zi is the position of the selected boundary. We number the layers and the 
boundaries such that the i-th layer is located at Zi-i < z < Zi. Then, the layers in the half-space z > Zi can be 
considered as passive (no radiation is comming out of them). In the circuit theory terms these layers constitute a 
load for the other, active, part of the structure located at z < z^, and can be equivalently represented by an input 
impedance, which, for a chain connection of four-pole networks, is given by the following recursive formula 

r^i+l yi+1 

^in+ - ^11 + 1 . + 2 ' ^^^-l 

^22 + ^in+ 
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where Z^^^ is the input impedance of all the layers behind the (z + l)-th layer. The recursion is terminated with 
the input impedance of the last layer which extends up to z = +00 (it can be, for instance, free space behind the 
structure), i.e., with the wave impedance of the last layer. Substituting ([T2|) into (|22|) for the two main polarizations 
in uniaxial layers we obtain 



■t+i 



where, for brevity, Zi+i — Z™'"""^ and /3. 



Z 



+ jZi+itan(^i+idi+i) 



i+l 



Z,- 



-iZ^+2tan(A+id 



(23) 



oTE.TM ^ 

. Eq. 



23)) is well-known in the theory of transmission lines. 
On the other hand, we may apply Thevenin's theorem to the active layers located at z < z^. Doing this, one first 
finds the internal impedance of Thevenin's equivalent circuit: 



ZL 



or, after substituting p^ . 



Z. 



22 



■2^12^21 



(24) 



Z.. 



Zl^} + jZt tan(/3i(ii) 

'z.+jZl-HMMrY 



(25) 



where Z^~_ is the internal impedance of the rest of the layers located at z < Zi-i, and the recursion terminates at 
the layer which extends down to z = —00. Note that because here we consider reciprocal structures, the Thevenin 
impedance ^f„_ equals the input impedance of all the layers located at z < z^ as seen by a wave incident from the 
half-space z > z^. 

Next, the equivalent voltage generator in Thevenin's theorem (recall that the EMF of this generator is the same as 
the output voltage of the network under the open circuit condition) can be found recursively as: 



^21 



^11 



Zl 



(26) 



where fg~^ is the equivalent EMF of all the sources located at z < Zi_i. This EMF is defined at the boundary 
z = Zi_i. Taking into account relations ([T5|) ~ (I^(I1) and the fact that the EMFs corresponding to distinct layers are not 
correlated, the mean square amplitude of fluctuations of £1. can be expressed after some algebra as 



^21 



z\i + zi 



2Re 



where 



'^21 



z: 



21 



Zl„ 



z: 



21 



Aw 



2R\^Q{^,T,)—^F^{^et'Y, (27) 

TT 



^li + Zl 



\Z, cos(/3,dO -f sin(;3,d,)P ' 



(28) 



and 



i?Jh=Re(Z?„_)~F*Re(Z-:i) 



(29) 



This important result shows that the effect of thermal fluctuations within the i-th material layer under the temperature 
Ti is fully equivalent to the effect of fluctuations in a resistance iJJjj placed under the same temperature. 

In other words, Eq. (f29| manifests that the input resistance Re(Z?j_) of a stack of layers can be split into two 
addends: Re(Z?jj_) — -I- F^Vi,e(Z\~^) . When considered together, these addends represent the total loss in the 
stack. However, when the thermal fluctuations are of concern, Eq. allows us to separate explicitely a part of the 
input resistance that appears in the Nyquist formula as being under physical temperature of the i-th layer. Thus, the 
noise produced within the z-th layer is associated with Rl^^. The other addend, F*Re(Z*^J'^), is due to the loss in the 
layers located below the z-th layer, and the thermal noise associated with it is understood as the noise received by the 
i-th layer from the background. 
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FIG. 4: Thevenin's equivalent network of thermal noise in layered media. 

Similar concepts exist, for example, in the antenna theory where the thermal noise of an antenna is represented as 
a sum of the noise generated locally by the ohmic loss in the antenna (analogous to RliJ and the noise received from 
the environment. The first addend in this case is proportional to the antenna loss resistance (which vanishes for an 
antenna made of a perfect conductor) and the second term is proportional to the radiation resistance of the antenna. 

From Eq. Re(Z?-i) = Rl^^ + F'-^ReiZi'^), therefore, we may as well write 

Rl^ = Re{Zl_) - F^Rl^' - F^F^-^Rl^' (30) 



where the series terminates at the layer (with the index i — M) that extends to z = — oo, for which RIy^^'^ 



Re (Zl^i- ^) = Re i^zj]^^^^ . We may analogously expand the last addend in (|27| which corresponds to the effect of 
fluctuations in the layers located at z < zi-i. Doing so we obtain 



(f^)2 = 2Rii,e{u;,T,)— + 2F'Ri^'e{u;,T,^,)— + 2F'F'-^ R',^^Q{u:,T,^2)— + (31) 

TT TT TT 



Thus, we conclude that the effect of thermal fluctuations in all layers located at z < is the same as in a chain 
of resistors with the values = R^^, R]^"" = F^R\^^, _R*^^ = F"" F'-''^ R\'^'^ , etc., kept under the temperatures Tj, 
Ti-i, Ti_2, etc. This result is analogous to the known formula for the thermal noise in cascaded amplifiers, in which 
case the quantities are called the noise factors. 

The corresponding equivalent circuit is shown in Fig. 21 in which we split Thevenin's internal impedance into a 
reactive part X^cn = Ini(^in_) and a resistive part Re(Z?j_) = X^^cff"- Respectively, Thevenin's EMF splits into a 



series of uncorrected fluctuating EMFs: f g = E ^off"' ^^i^^ (^off")^ = (2/7r)e(a;, ri_„)i?^g"Aw, n = 0, 1, 2, . . . , M, 

n 

representing the effect of thermal fluctuations in the layers located at z < . The rest of the structure at z > z^ is 
modeled by an effective load impedance .^load = ^l^^. 

The radiative heat flux from any layer located at z < z^ into the half-space z > Zi can now be trivially calculated 
based on this equivalent circuit. Namely, the power spectral density associated with the plane waves with a given 
transverse wave vector kt — {kx,ky) and a given angular frequency u, originated in the slab with the index i — n, can 
be expressed as 



^c^.kt - TvTTTi i^Re(Zioad) - -— , y,+ l|2 ^^i4n+j • (32) 

Respectively, for the total radiative heat flux (associated with waves of a selected polarization: TE or TM) into the 
half-space z > z^ we have 



oo 



"00 

The flux crossing the the same boundary in the opposite direction is found by reversing the roles of the active and 
passive layers. 

IV. PARTICULAR CASE I: BLACK BODY RADIATION 



Let us apply this equivalent circuit theory to calculate the power radiated by a black body per unit of frequency 
and surface area. We assume that a very thick black body lies in the lower half-space z < and is under the constant 
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temperature T. The upper half-space z > is empty. We are interested in the thermal radiation into this half-space 
from the black body surface at z = 0. 

The equivalent circuit for this system is composed of a single resistance R^^^ = Re^Zj^j^l^, where z\^_ is the input 
impedance of the half-space z < occupied by the black body, a corresponding fluctuating EMF , a reactance 
Xgcn — Iml ) , and a load Z\oad = ^in+ — ' ■> which is the input impedance of the open half-space at z; > 0. 

By definition, the black body absorbs all incoming radiation independently of the frequency or the angle of incidence. 
Thus, electromagnetically, there are no reflections from such a body which means that it is perfectly impedance- 
matched to the free space. Therefore, z[^2- = with 



^0 — 



Vmo/£o 



/7TM _ 
^0 — 



(34) 



and, using Eq. p2p . we may write the power spectral density associated with the radiative heat flux into the open 
half-space as 



2e(cj,T) 



Re 



2Z, 



TE,TM 



0, g > fco 



(35) 



The result is the same for both polarizations. From here, the total power emitted from the black body by both 
polarizations per unit of its surface, per unit of frequency is 



— - 2 X 



(2^) 



Q{lo,T) 

47^3 



dk^ dky 



(36) 



fe^+fe2<fc2 

The same result can be, of course, obtained from Planck's expression for spectral radiance of a black body which 
reads 



■e(w,r). 



(37) 



The spectral radiance is defined as the power emitted from the black body surface per unit projected area of the 
emitting surface, per unit solid angle, per frequency: 



APr 



rad 



Aj_AnAu} 

Hence, because in spherical coordinates the projected area Aj_ — AqCO'&B^ we find 



dS_ 
(ku 



2iT /.tt/2 



B^{T) cos 9 sine d0 d(p 



'0 JO 

which is the same as the result predicted by the equivalent circuit model 



A'k'^c 



re(w,r). 



(38) 



(39) 



V. PARTICULAR CASE II: GENERALIZED POLDER- VAN HOVE FORMULA 



In this section we derive a generalization of the Polder- van Hove formula applicable to layered uniaxial magneto- 
dielectrics, using the equivalent circuit model developed in Section iHll The geometry of the structure is the same as in 
Fig.[T](a). We are interested in the radiative thermal transfer between the media which occupy the half-spaces z < 
and z > d (the media with indices 1 and 3). These half-spaces are kept under temperatures Ti and T3, respectively. 
The region < z < d (the region 2) may be filled with another uniaxial medium kept under temperature T2, or may 
be left empty (which is the case of a vacuum gap). 

In order to find the radiative heat flux from medium 1 into medium 3 we split the structure at the plane z = d, and 
consider the layers located at z < d as active layers. The half-space z > d plays the role of a load. The equivalent 
circuit of such a structure can be represented as in Fig. HI with a pair of resistors ii^g^ = F'^^^r[1^ = F^'^^Re{Zi) and 
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■^eff ~ ■'^^(■^in-) ^ ^cff ' ^ P^i^ °f corresponding fluctuating EMFs and , a reactance Xgcn = I™(^^in-) j 
and a complex load Zioad = ■^3- 

(2) 

One may verify that in the vacuum gap case R^g = 0, which is a consequence of the fact that there is no dissipation 
in the gap. Evidently, the same conclusion holds when the gap is filled with a lossless medium. Nevertheless, 
the following derivation is general enough to be applicable to both medium-filled or vacuum gaps, with or without 
dissipation. 

We start with calculating the noise factor Because Z^^^l = Zi, we obtain from ()28p : 



K2) _ 



|Z2COs(/32d)+j^isin(/32d)P |1 - ri2e-2j&rf|2 ' 
where ri2 is defined in Section |lll Thus, from (|32|) . 



Next, Thevenin's internal impedance Z^^J_ is found from (1^ 



(40) 



(2) ^ Zi+jZ2tan(/32d) ^ 1 + ri2e-^^fe'^ 

'Z2+jZitan(/32d) ^l-ri2e-2j732d' ^^^^ 



from which the total impedance of the network reads 



^(2) ^ ^ _ 2Z2 i-ri2r32e"^^fe'^ 



where r32 is defined in Section |TI1 Substituting (|40l) and into (HIT) we obtain 

^ 6(0., Ti) |l-ri2ni-r32p|e-^fe'^p Re(Zi)Rc(Z3) 

2^ |l_ri2r32e-2j/32rf|2 |Z2|2 ' ^ > 

Respectively, the total radiative heat flux associated either with TE or TM polarized waves originated from medium 
1 and absorbed in medium 3 is 

CXD CXD 

^1^3 = j^JdioJ qdqe{u,n)N{cu,q), (45) 



where 

^ '^^ |l-ri2r32e-2j/32d|2 4|Z2|2 ^ ^ 

In the vacuum gap case, the wave impedance Z2 and the propagation factor /32 are purely real (imaginary) for the 
propagating waves (evanescent waves) in the gap. Therefore, because Zi = ^2(1 + ri2)/ (1 — ^12), * = 1, 3, we have 

1 — IF- P 

I ^2 1 1 -I „ / 1^ ' (propagating waves) , 

MZ.) = <! (47) 

^\Z2\tz (evanescent waves), 

|1 - 1j2| 

which results in the Polder- van Hove formulas when substituted into (j45l) and (146^ . However, the more general result 
represented by Eqs. (j45p - (|46|) holds for arbitrary uniaxial magneto-dielectric media filling the gap. It is easy to verify 
that Eq. pS)) can be as well written in form ©. 

Moreover, in general, when the gap is filled with a lossy medium and T2 7^ one also has to take into account 
the radiative heat flux to medium 3 that is originated in the gap (i.e., in medium 2). The power spectral density 
associated with it is, from Eq. ([32]), 



2 Q{u,T2)R^^^ e(^,T2)|l-f 12^11 -r32p|flfff^Re(Z3) 
^uj,kt - ~~l2\ r^'S Z3J - ^— — — , (48J 

^\z\:l + z-iV 27r |i-ri2r32|2 \Z2Y 
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where the reflection coefRcient defined at the plane z = d. One may note that (1481) has the 

same form as (|44|) with d = and one of the wave impedances replaced by the effective resistance ■ This is the 

consequence of the fact that when the thickness of the middle layer increases, — > 0, Rl^f^ Re{Z2), and (gll) 
reduces to the Polder-van Hove's result for two media in direct contact. 



VI. ANTENNA THEORY AND CIRCUIT THEORY CONCEPTS APPLIED TO RADIATIVE HEAT 

TRANSFER 

In this section, in order to better understand the circuit model developed in this work, we establish a connection 
between our model and the classical theory of noise in receiving antennas. We consider the radiative heat transfer 
example from Section |V] and employ an analogy between the heat-receiving half-space (medium 3) and a loaded 
receiving antenna. More exactly, in this analogy the unit area of the interface between the media 2 and 3 is treated as 
an aperture antenna that receives the power of thermal radiation from the half-space z < d and delivers it to medium 
3 which is understood as the antenna load. In what follows, we assume that medium 2 is lossless, therefore, all the 
radiative heat delivered to medium 3 is generated in medium 1. 

The equivalent circuit of this problem is the same as the one discussed in Section |V] with i?^^ = 0. Thus, in the 
antenna analogy there is one noise source with the internal impedance ZIJ_ = Za = Ra + i^A, where Xj^ — X^cn is 
the antenna reactance and i?A — -F'^^^Re(Zi) is the radiation resistance of the antenna. Such analogy signifies that 
the effective radiation resistance of an aperture equals the real part of the input impedance of the half-space seen 
from the aperture. The antenna load is represented in this analogy by the impedance Zioad = -^3 which is the wave 
impedance of medium 3. 

As is known from the theory of noise in lossless antennas^ (an aperture by itself has no loss), the mean-square EMF 
of the thermal noise of a directive antenna (e.g., a radio telescope) is given by the Nyquist formula ([T]), in which 
one inserts the radiation resistance of the antenna [as R{y)] and the effective temperature of the area of the sky to 
which the antenna is directed (which is medium 1 under temperature T = Ti in our analogy). Applying this to the 
antenna equivalent circuit, we may write for the noise power at the antenna load: 

^out = ^^Re(Z3) = -^P^Re{Z,). (49) 

|za + z-iY TT |Za + z-iY 

For our example of a lossless medium 2, i?A = Rc {z'^^^ , which brings us to the same result as the more general 

cascade-circuit model (HIT) , i.e.. Pout = Pli^- Let us also note that if there would be no separating layer, then Z^ 
would be simply equal to the wave impedance of medium 1. With the separation layer in place, antenna model calcu- 
lation is equivalent to calculation of input impedance of a transmission-line section loaded with a known impedance 

[Eq. m\- ^ ^ ^ 

It is worth noting that the tight connection between the model of the present paper and the theories of noise in 
antennas and cascaded electric networks allows one to better understand optimal conditions for radiative heat transfer 
through composite layers and, consequently, to design these material structures aiming for desired and optimized 
performance. In particular, from (1491) and (|4ip we see that the problem of maximizing radiative heat transfer for 
given layer temperatures reduces to an equivalent problem of matching a generator to a load. Let us discuss this issue 
assuming for simplicity, that the media 1 and 3 are the same, i.e., Z\ = Zj,. At the first glance, it appears that the 
optimal heat transfer is ensured if we simply connect the two equivalent media together (or fill the gap with the same 
medium as those 1 and 3). However, this is true only if the wave impedance is real. For complex Z\ and ^3 (which 
is realistic even for propagating modes in view of losses in the media), the best radiative heat transfer corresponds 
to the conjugate impedance matching Z\ = when the negative reactance of one of two media is compensated by 
the positive reactance of the other one. Then, the spatial spectrum of the heat transfer function in accordance to ^ 
turns to -/V = 1/4 whereas for the direct contact of two equivalent media with wave impedance Z = R-\- jX we have 
N = R? I\2ZY ■ This tells us that the radiative heat transfer through a properly filled gap can be in principle made 
larger than that through the direct contact of two equivalent media. 

In Ref. Isl it was proposed to insert a slab of the so-called indefinite medium between media 1 and 3 [Fig. [T] (b)]. 
Indefinite media (also called hyperbolic metamaterialsiS) are uniaxial dielectrics characterized by the permittivity 
tensor that has opposite signs of the longitudinal and transverse components. Such filling allows for enhancement 
of the heat transfer by increasing the noise factor F^^^ . This factor increases compared to the vacuum gap because 
indefinite media support propagation of spatial harmonics with high transverse wavenumbers, which would otherwise 
be evanescent in the gap [notice the exponential factor in Eq. (001)] ■ For such waves, F'^' dramatically increases at the 
spatial frequencies that correspond to the minima of the denominator of Eq. (I40[) . However, it is not straightforward 
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to ensure proper impedance match in such structures. In this view, the structure suggested and studied in Ref. 
is not fuUy optimal, although it still demonstrates that a specifically crafted filling may dramatically enhance the 
radiative heat transfer through the gap. 

The transmission line analogy, however, suggest an immediate possibility to circumvent the problem of impedance 
conjugate match. The key is to make the wave impedances of all three media real (at least, approximately) and 
equal in all layers. This can be realized using a wire medium in region 2 which extends inside regions 1 and 3. 
Indeed, it is known that wire media support propagating TEM modes with high spatial frequencies (transverse wave 
numbers) q, including q > ko and limited only by the period of the wire array (see, e.g., in Ref. 23i). Wires should be 
good conductors at the frequency of interest, and the background media nonconducting. If the wires extend over all 
three regions and none of them has a conducting background, then these TEM modes will exist and have real wave 
impedance everywhere in the system, in principle allowing the conjugate match of the load to the heat source for a 
wide range of spatial frequencies q. 

In the present paper we, however, do not present the studies of the heat transfer optimization keeping these results 
for our next publications. Instead, in the next section we report some numerical results illustrating the applicability 
of our model for solving practical problems. 



VII. NUMERICAL EXAMPLE: RADIATIVE THERMAL TRANSFER THROUGH A 

NANOSTRUCTURED LAYER 



In order to demonstrate the applicability of the developed theory to real- world problems, we calculate in this section 
the spectral density of the radiative heat flux absorbed in medium 3 of the structure depicted in Fig.[l](b): 

oo 


We compare the heat transferred across the vacuum gap with that transferred across the gap filled with either 
normally oriented metal-state single-wall carbon nanotubes (CNT) or with similarly oriented golden (Au) nanowires. 
Calculation for the array of CNT is done in order to validate the present model using the exact simulations of Ref. 
Calculations for metal nanowires are done in order to confirm or decline the effect of giant enhancement of radiative 
heat transfer in the near infrared (IR) range due to the presence of nanowires. This effect was predicted but not 
studied in Ref. H. In accordance with the theory presented above, all calculations are done for homogenized media. 
The homogenization models for j-th medium results in explicit formulas for e,^ and (whereas = ^| = 1). 

The mid-IR homogenization model for the arrayof CNT was described and validated in Ref. [2J. The parameters 
of the array of CNT correspond to those of Ref. [a [see also in Fig. [T] (b)] . The homogenization model for aligned 
metal nanowires representing their array as a layer of an indefinite material was described in Ref. [25l . It is applicable 
to the visible and near-IR ranges where it offers a rather high accuracy for optically dense arrays of rather thin wires. 
Practically, for the band A =1-2 /im one needs the period below 300-600 nm and the wire diameter of the order of 
the skin depth in the metal or smalleri^ For Au in this wavelength range it means that the thickness of the wires 
should not exceed 20-50 nm. 

In both cases (CNT in the mid-IR range and nanowires in the near IR) the presence of the metamaterial enhances 
the radiative heat transfer. This effect results from the conversion of TM-polarized evanescent waves into propagating 
ones in the effective indefinite material filling the gap.— i2ii2^ Respectively, in this section we analyze only that part of 
the radiative heat which is transferred by TM-polarized waves. 

In this numerical example we neglect the contribution of thermal sources located in medium 2 (i.e., in CNT and 
nanowires). The impact of thermal production and absorption in medium 2 will be evaluated in our next paper, where 
we will also consider thermo-photovoltaic applications of our present model. As is clear from Fig.[T](b), medium 2 in 
the gap is a triple-layered structure, therefore, there are in total five material layers in the whole structure. However, 
since in this example we neglect the thermal processes in the gap, we may replace the layers in the gap by a single 
equivalent four-pole network. Its Z-matrix, Zmn, is obtained in a standard manner from its transfer matrix^^i with 
the latter being a product of transfer matrices of effectively homogeneous anisotropic layers with thicknesses hi , h2 
and /i3. Next, the power spectral density P^~^^ in Eq. ([5(1)) is calculated from Eq. (PT|) where the noise factor is 

= Z^l^ ^ I Z[f + Zi ^ [Eq. dMl)], and the internal impedance is Z^^l = Z^l^ - z[lh^f / (^z[f + Zi j [Eq. (HH)]. 

In Fig. [5] (a) we depict the energy transfer coefficient |rp introduced in Ref. i [Eq. (12) of Ref. @ which differs 
from our heat transfer spatial spectrum N{uj,q) defined by Eq. (fTU)) . by the factor 4|Z3p/i?ii?3. In the present case 
media 1 and 3 are equivalent (heavily doped silicon), i.e., Zi = Z3, Ri = R3. All parameters in the calculation 
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— Vacuum gap, exact 
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No nanowires 




(a) 



(b) 



FIG. 5: (Color online) (a) - Coefficient |rp — 4|Zi/_RipA^ at A = 7.5 versus normalized spatial frequency qa/7v for the gap 
filled with free space and for that filled with an array of CNT. Dashed lines - exact (beyond homogenization) simulations from 
Ref. IB- Solid lines - analytical calculations in accordance to the present model, (b) - Spatial spectrum A'^ of the radiative heat 
transfer function versus dimensionless spatial frequency q/k at A = 1.5 fj,m. Calculations are done for four cases - nanowires 
are only in medium 2, nanowires are in both media 1 and 2, nanowires are in medium 1 only and nanowires are absent. 
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FIG. 6: (Color online) The gain in the spectral density of the transferred radiative heat flux due to the presence of nanowires 
(two cases of their arrangement), (a) - Gap of thickness d = 2 ^m. (b) - Gap of thickness d = 0.5 fim. 



illustrated by Fig. [5] (a) correspond to those from Ref. S The value |rp is calculated at wavelength A — 7.5 /im 
as a function of the normalized spatial frequency qa/n, where a = 20 nm is the period of the CNT array in the 
domains hi = d/3 and — d/3, where the gap thickness is d = 1 /im. In the domain /12 = d/S the array period 
is equal 02 = \/20 nm. Heavily doped silicon supports so-called surface plasmon-polariton (SPP) waves generated 
on the surfaces of the Si half-spaces at {qa/ir) = 0.01 in the case when the gap is empty. The value q = O.OItt/o 
nearly corresponds to g = l.Sfco, where fco is the free-space wave number. The manifestation of this SPP is the 
maximum of the corresponding curve in Fig. [S] (a). Dashed curves in Fig. [5] (a) correspond to exact simulations of 
Ref. m which take into account the microstructure of the material layer in the gap. In Fig. [S] (a) we show only the 
region of spatial frequencies (qa/n) < 0.08 in which the difference between the exact and homogenized models of the 
CNT array is negligibly small (see Ref. For both empty gap and gap filled with CNT the agreement between 
our model and the exact simulations is excellent. Local maxima of the transmittance spatial spectrum for the gap 
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filled with CNT correspond to the thickness resonances of spatial harmonics [in presence of CNT the whole region 
{qa/n) < 0.08 corresponds to propagating plane waves, though the inequality q > holds for {qa/n) > 5.6 • lO^'^]. 
Fine agreement between our circuit model and exact simulations pertains at other wavelengths besides the ones 
mentioned in Figs. [5] (a,b), because with our circuit model we also have reproduced the frequency dependence of the 
gain in the heat transfer G{ui) = s'\^^'^\uj) / sf^ (tu) , calculated in Ref. 's'. Here sf^ corresponds to the vacuum gap 

and s^^^"^^ corresponds to the gap filled with CNT. 

In Fig.[S](b) we present the dependence N{q, A*) for the case when the gap d = 2 fj,m is filled with golden nanowires. 
The complex permittivity of gold in the range A = 1 ... 2 /j,m was taken from Ref. Regions 1 and 3 in this case 
are filled with doped germanium used in real thermo-photovoltaic systems whose complex permittivity was taken 
from Ref. [13. Function N{q, A*) was calculated for nanowires with volume fraction p — 0.2 in the domains hi and 
hs and p = 0.4 in the domain h2 (when hi = h2 = hs). Here A* = 1.5 /im has been chosen having in mind possible 
thermophotovoltaic applications (at this wavelength the doped Ge has nearly maximal photovoltaic spectral response). 
The result for N (thick solid curve) was compared with that for the empty gap d = 2 /.tm (thin dashed curve) at the 
same wavelength. 

It has to be mentioned that the structure shown in Fig. [T] (b) with free-standing metal nanowires is an abstraction 
— in a feasible structure nanowires are partially submerged into the host material. Therefore we have also calculated 
the function N(q, A*) for the case when Au nanowires are semi- infinite and have the same volume fraction p — 0.3 in 
media 1 and 2. This calculation [thick dashed curve in Fig.[5](b)] is done in order to understand how the extension of 
nanowires into medium 1 changes the radiative heat transfer. In this case the gap is uniformly filled, i.e., in the the 
structure shown in Fig. [T] (b), h2 — h^ = and hi = d. In this case the nanowres touch the surface of medium 3 and 
the thermal transfer by the direct thermal conductance may be of significance, in addition to the radiative one. This 
effect is not considered in the present paper. Additionally, we have studied the case when nanowires with p — 0.3 are 
located only in medium 1 and the gap is empty [thin solid curve in Fig. [S] (b)]. 

We can see in Fig. [5] (b) that the integral increase of N{q) with respect to the empty gap is very significant for 
both cases when the nanowires fill in the gap. Function N(q) in the case of the vacuum gap has two local maxima in 
the region q < ko resulting from Fabry-Perot resonances. Because the structure is not fully impedance-matched these 
maxima are much smaller than the achievable limit = 1/4 and N{q) vanishes fast at g > feg. Unlike the situation 
illustrated by Fig. [5] (a) the real part of the complex permittivity of Ge at A = 1.5 /xm is positive and SPP cannot 
be excited. Since the value of N{q) for the empty gap is so small, the gain granted by Au nanowires to the heat 
transfer between two half-spaces of Ge turns out to be larger than that offered by CNT to the heat transfer between 
two half-spaces of Si calculated in Ref. d. 

The presence of nanowires only in medium 1 turns out to be destructive for the amplitude of N. In this case 
medium 1 is an indefinite metamaterial, and the mismatch between medium 1 and free space increases. However, if 
nanowires are present in both media 1 and 2, noticeable values of N keep for q > ko. The case when nanowires are 
located only in medium 2 is the best one: it corresponds to the smallest mismatch between the media. 

The dependencies shown in Fig. [5] (b) are typical for every A in the band of the photovoltaic operation of Ge 
(A = 1 — 2 fim). As a result, due to the presence of nanowires the heat transfer gain G = Sia'^V'^ia'' almost uniform 
over a wide range of wavelengths. Here s\^^^ corresponds to the case of nanowires in the gap and s^i^ corresponds to 
the vacuum gap. In Fig. [6] (a) we present the gain in dB, i.e., 10 log^g G calculated for the case d = 2 iim. The huge 
gain keeps for the interval of values d = 0.5 ... 5 fim. In Fig. [6] (b) we show the same gains for the case d = 0.5 /im. 
We can conclude that the presence of Au nanowires in the micrometer or submicron gap can increase the near-IR 
energy transfer across the gap by 3 orders of magnitude. This result confirms the expectations of Ref. 



VIII. CONCLUSIONS 



In this work we have formulated an equivalent circuit theory of the radiative heat transfer in uniaxial stratified 
magneto-dielectric media. We have proven that the effect of thermal-electromagnetic fluctuations in such structures 
can be fully determined without an explicit knowledge of the microstructure of the layers, as well as without a need 
to employ any calculations based on distributed fluctuating currents. Instead, the only physical characteristic on 
which we base our theory is the effective input impedance of a stack of layers, which can be obtained for any spatial 
harmonic of the field (including both propagating and evanescent waves) using the methods of VTLT. We have shown 
that such impedance representation, while being in full agreement with sophisticated full-wave methods known from 
the literature, results in simple formulas analogous to Nyquist theory-based formulas for thermal noise in cascaded 
electric circuits (for example, cascaded amplifiers). Therefore, with this model some important concepts from the 
theory of electric networks (conjugate-impedance match, optimal filtering, etc.) can be imported into the field of 
radiative thermal transfer in multilayered structures. 
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From the point of view of practical implementations, the developed equivalent circuit approach offers significant 
simplifications as compared to the known theories of radiative heat transfer based on distributed fluctuating currents. 
Without any modifications, our method can be used in heat transfer studies in uniaxial anisotropic media that include 
micron and (or) submicron-thick layers. Moreover, our model is applicable to radiative heat transfer in composite 
or nanostructured layers (if these layers are effectively homogeneous for spatial harmonics of the electromagnetic 
field which transfer the radiative heat), and is readily generalizable to stratified bi-anisotropic and spatially dispersive 
materials. Therefore, we hope that our work may significantly enlarge the scope of the radiative heat transfer research 
in composites, especially in nanostructured metamaterials. We believe that this may lead to new opportunities in the 
design of efficient thermal energy harvesting devices, like thermophotovoltaic converters and such. 



Appendix: Nyquist formula for a reciprocal anisotropic and lossy magneto-dielectric slab 

We consider a uniaxial magneto-dielectric slab described by the macroscopic Maxwell equations for the time- 
harmonic fields 

V X E = -jcjp^ • H - J", V X H = jcjfa-E + J", (51) 

with the absolute permittivity and permeability dyadics of the form Sa = £o(£^-^t + e"zoZo) and /la — fJ-oil^^It + 
/i'lzozo). We assume that the material of the slab is lossy, therefore, by the fiuctuation-dissipation theorem (at 
non-zero temperature) there appear fluctuating external currents J*^ and J'" in the slab. The explicit form of these 
currents is not important at this stage. As in the main text, here we use the convention in which the time-harmonic 
quantities are understood as root mean square (rms) values. 

Let the fields E', H' be an arbitrary solution of the Maxwell equations ([51]) with J*^ = J'" = within the slab. Then, 
considering the two systems of Maxwell equations with non-zero sources and with vanishing sources, respectively, we 
can form the Lorentz lemma 

V • (E X H' - E' X H) = E' • J'^ -H' • J". (52) 
Integrating it over the volume V of the slab, we obtain the reciprocity relation 

f ni • (E X H' -E' X H)|^^ d5+ / n2 • (E x H' - E' x H)|g^ rfS* = / (E' • J-^ - H' • J") dV^, (53) 

J Si 'J S2 V 

where Si^2 are at the two interfaces of the slab, and ni.2 are the outer unit normals to these surfaces, respectively. 

One may select any solution of the uniform Maxwell equations within the slab for the fields E', H'. For us it is 
convenient to use the one that has the form 

E'(r) = E'_k,(2)e^'^*-, H'(r) = H'_k, (^)e^■'''•^ (54) 

where the z-axis is orthogonal to the slab and the real vector kj lies in the plane of the slab (the xy plane). Physically, 
such a form corresponds to a superposition of plane waves with the same transverse wavenumber: — kf. In Eq. (j54l) . 
E'_i^^ and H'_jj.^ define the field solution profile within the slab as a function of z. 

Substituting (|54|) into the reciprocity relation (|53|) . we obtain (the first slab interface is at 2; = zi and the second 
one is at z = Z2) 

ni • (Ek, X H'_k, - E'_k, X HkJl^^^^ + ■ (Ek, x H'_k, - E'_k, x HkJ|^^^^ = 

''(E'_k, -J^.-H'^k, -Jirjd^, (55) 

1 

where we have decomposed the fiuctuating currents J*^-™ and the fields E, H into plane waves using the Fourier 
transform defined as 

F(r) = (^// Fk,(z)e-^'^''-d2k„ Fk,(z) = ^// F(r)e^''-"- d^r, (56) 

where Aq is the unit area in the a;?/-plane, and F can be any of the fields or currents. 

Eq. (|55p is the reciprocity relation for the wave components characterized with a fixed transverse wavenumber. In 
order to simplify further writing we will use the notation Fi 2 = F±kt (2:1,2) with F being any of the fields or currents. 
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Then, noticing that only transverse components of the fields play any role on the left-hand side of ([55]) we rewrite it 
as 



ni • (El, X h;, - E'l, X Hi, 



n2 • (E2, X H^, - E^, X H2J = 

E'l, • (m X Hu) - El, • (m X H'l,) + E^, • (n2 x H2,) 



E 



2t • 



(n2xH^,). (57) 



Let us remind that the quantities E'l 2^ and H'l 2^ have the meaning of the transverse components of the electric 
and magnetic fields at the interfaces of a source-free magneto-dielectric slab. Therefore, as follows from the vector 
transmission line theory (VTLT) for such slabs, these components are related by the impedance matrix of the slab 



E'l, \ I Zxx Z12 
2t J \ Z21 Z22 



E 



ni 
n2 



H'l, 



H' 



2t 



(58) 



The components of this matrix are dyadics that are even in k,: Zmn{—^t) = ZmnO^t)- Also, due to the symmetry 

= =T z 

12 = ^21, and Z = . 



and the reciprocity, Zn = Z22, Z 
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Using ([55)1 on the left-hand side of Eq. ([57)1 we obtain 



^0 E'l, • (m X Hi,) - El, • (m X H'l,) 



E^, • (n2 X H2,) 



I'l- Zii-Ii 



E 



2t • 



^12 



(n2 X H^,) 



Vi 



I'. 



^21 ■ Ii — V2 



(59) 



where we have introduced the vector currents I'l 2 = -\/^ni 2 x H'l 2^^ Ii,2 = \/~^t^i.2 x Hi^2, and the vector voltages 
Vi_2 = a/Aq Ei_2,, and used the symmetry properties of the impedance dyadics. 

Let us now work on the right-hand side of ([55]). At a fixed k, the Maxwell equations for the fields E'(r), H'(r) 
reduce to a system of first-order linear differential equations for the vector functions E'_jj.^(z) and '^^'_-^^^{z). From the 
uniqueness theorem it follows that these functions are univocally defined by boundary conditions imposed either on 
tangential electric or tangential magnetic field. Thus, we may consider two auxiliary boundary-value problems, the 
first one with the boundary conditions 



and the second one with 



m X HLk,(^i) - I'l/v^, n2 X HLk,(22) = 0, 



ni X HL. (zi) = 0, n2 X HL. (^2) = Iz/V^- 



(60) 



(61) 



The field equations are the same in these two problems. From linearity it follows that the superposition of the solutions 
of the two problems is the same as the fields E' ^.^(z) and H'_j^^(z) that appear in ([SQ)) . On the other hand, these 
problems physically correspond to the two cases of the magneto-dielectric slab backed with a magnetic wall (perfect 
magnetic conductor, PMC) at z = Z2 and at z = zi, respectively. 

Let us consider the problem with the boundary conditions (|60p . We may split the vector I'l into the components 
parallel and orthogonal to k,: 

T' _ r' ^* _L r' ^* ^ "1 

■^1 — ^l,TMTi — r+^l,TE — Tj — i • 



Thus, the component I[ rpj^ corresponds to TM-polarized field, and the component I[ corresponds to TE-polarized 
field. Because the wave equations in the slab also split into independent equations for the TM and TE waves, we may 
also write for the vector fields 



ELk,(^) = E'i^tm(^) + E'i,teW, HLk,(^) = H'i,tm(^) + H'i,te(^), 



(63) 



where the addends are the TM and TE solutions for the fields in the PMC-backed slab. The magnitudes of these 
solutions are proportional to I[ rpj^ and I[ rpg, respectively. 

Based on the above discussion we find for the right-hand side of ((55)) 



(eL 



Ul^^-3l)dz 



pZ2 PZ2 



A, 



■I'l 



ei,TMi 



ei,TE — Tj— j — , (64) 
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where 



ei,TM 



ei,TE 



TM -J zi 



Yi / (^l.TE ■ Jkt ~ • JjfJ dz. 

''l,TE -I zi 



In an analogous manner we consider the second case with a PMC bX z — Zi and find 
where 



/' 



^-kt • Jkt - il-kt ■ -iicjdz = —L^- |^e2,TM + e2,TE— 



e2,TM = / (E'a^TM • Jkf - H2,TM • ^kt)dz, 

^2.TM J z. 



(65) 



(66) 



(67) 



(68) 



62, TE 



^0 /"^^ 

77 / (E2,TE ■ Jkt ~ H2_xE ■ JlTt) 
-'2,TE -Jzi 

Therefore, combining these results together and using (|55p . ([57l) . and ([59l) we obtain 



I'l 



Zii • Ii + Z12 • I2 — Vi 



Z21 • Ii — V2 



,/ / kt X ni 

Il-(ei,TM^+ei,TE^^ 



I-) • I e2,TM7; — r 



Finally, because the vectors ![ and I2 are arbitrary, 

\ Z21 Z22 ) 



Vi 
V2 



ei 

62 



(69) 



62, TE- 



kt X n2 



(70) 



(71) 



where 61^2 — ei,2,TM(kt/|kt|) +ei,2,TE(kt x ni_2)/|kt|. These equations represent the equivalent vector circuit model 
of a magneto-dielectric slab with fluctuating sources. In this model, Ii,2 have the meaning of equivalent vector currents 
at the two ports of a linear four-pole network of dyadic impedances, and ei_2 are the equivalent vector EMFs acting 
at the two ports. 

Because the equivalent EMFs are expressed through the fluctuating currents, they are also fluctuating, stochastic 
quantities. As is readily seen from ((65)) - (|66)) and (l68t - ([69|) . the stochastic mean value of the fluctuating EMFs is 
zero: ei/2 — 0, because 3'^^ = 0. However, the mean-square values of the fluctuating EMFs, as well as their mutual 
correlations are in general different from zero and can be calculated as follows: 



«,pe/3,9) 



TTT^ f\^'l,, - Jkt* - H';, • dz r\wp., ■ H'^,, . J^) d.' = 

^ a,p^ f},qJzi J zi 

jpAr- r r (^'".p -^V- ^'Ip ■ JkDL • J£t - HV,, • dzdz' 

^ a,p^ P,q Jzi Jzi 

a2 r I-Z2 pZ2 

° / / ^%{z)-3i:{z)3lSz')-^'p^,{z')dzdz'+ 

^ a,p^ P,q Uzi J zi 

r rn'l^iz)-3^;iz)3^^iz')-Il'p^^iz')dzdz' 



, (72) 



where a,/3 = 1,2 and p,q = TE, TM. There are no cross terms in the last integral of Eq. (|72l) because the electric 
and magnetic fluctuations are statistically independent: the dyadic Jk*Jirt ^^^^^ ^^^^ Jk^JiTf ^ ^- 
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From the fluctuation-dissipation theorem, for the fluctuating currents composed of harmonics within a narrow 
interval around a given frequency w, 



-(z)J£^(z') = -i-jc. (sa - el) 5{z - z')Q{lo,T)^u 



(73) 
(74) 



The dimensionality factor 1/Aq appears in ([75)) -(fM l) because of the form of transformation ([55]). Note also that ([75)) - 
([74| are written for the rms amplitudes of the fluctuating currents. 

Substituting ([75)) - (rMl) into (|72l) and evaluating the integrals over z' we find that 



jujAoQ{Lo,T)Auj 



E' 



'a,p ■ - el) • E'^^q + H'^ ^ 



a A^a 



Ma - A^aJ ■ H'/3,g 



However, from the well-known differential lemma 



V • (El X H* + E* X Hi) = -jco [e; ■ (ea - el) • Ei + H* • (/I, - Jll) ■ Hj 
which holds for arbitrary source-free electromagnetic fields Ei 2(r), Hi 2(r) within the slab, it follows that 



(75) 



(76) 



E' 



'*a,p ■ (ea - e^) • E'^^g -I- H'* p • (^^l^ - ^^) • li'p,q 

1 



dz 



n2 • (E'^,, X H'lp)|^ + 112 • (E'lp X H'^,,)|,^, 



(77) 



from which we see that ii p q, the integral (177)) vanishes due to the orthogonality of the TE and TM polarizations. 
Next, when p^q and a = ^ = 1 we obtain from ([77]), (inD])-(IlT]), and ((55)) : 



E'l p • ^ea - cl) • E'l^p + H'l p • 



• (Ma-Aia) -H'l 



az = f — 



Re(Zfi). 



(78) 



An analogous result is obtained ior a — /3 — 2. On the other hand, when a — 1 and /3 = 2 we obtain 

z = :^^2ii(zf2 + zfi*). 

Combining all these results together and using the reciprocity property of the Z-parameters we find from ()75p that 



J E'l p • - e^) • E'2^p + H'l p • (^/i^ - fil^ ■ H'2,j 



(79) 



(e*,pe^,p) = -RefZ^Je(w,r)Aa;, 



(80) 



where a,/3 — 1,2 and p = TE, TM. Eq. (jSO)) is the generalized Nyquist formula for the thermal-electromagnetic noise 
in a uniaxial magneto-dielectric layer. 
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